**************************************************************
* ELECTRICITY SPILOVERS PROJECT
*  MAIN REGRESSIONS
**************************************************************
set more off
frame copy data_main main_regs
frame change main_regs
keep if post==1

********************************************************************************
*TABLE 2: PANEL A - ELECTRICITY INTENT TO TREAT EFFECTS
********************************************************************************
eststo clear
*Column 1: Full Diff + weather + date + HOD FE
	eststo: reghdfe euse ws temp_60 temp_65 temp_70 temp_75 temp_80  precip_1, ///
			absorb(hour date) vce(cluster HHID)
	sum euse if e(sample) & ws==0
	scalar euse1 = r(mean) 
		estadd scalar m_cont = scalar(euse1)
		estadd local weath "Yes"
		estadd local hod "Yes"	
		estadd local date "Yes"
		estadd local hours "All"
		estadd local sam "05/2015-05/2016"		
	est store A
	
*Column 2: Full Diff + weather + date + HOD FE + pre-treatment use
	eststo: reghdfe euse ws elec_use_summer elec_use_annual elec_use_winter temp_60 temp_65 temp_70 temp_75 temp_80 temp_85 precip_1, ///
			absorb(hour date) vce(cluster HHID)
	sum euse if e(sample) & ws==0
	scalar euse2 = r(mean)
		estadd scalar m_cont = scalar(euse2)
		estadd local weath "Yes"
		estadd local hod "Yes"	
		estadd local date "Yes"
		estadd local hours "All"
		estadd local sam "05/2015-05/2016"	
	est store B
		
*Column 3: Summer 2015 Diff weather + date + HOD FE 
	eststo: reghdfe euse ws temp_60 temp_65 temp_70 temp_75 temp_80 temp_85 precip_1 if ym<=667, ///
			absorb(hour date) vce(cluster HHID)
	sum euse if e(sample) & ws==0
	scalar euse3 = r(mean)
		estadd scalar m_cont = scalar(euse3)
		estadd local weath "Yes"
		estadd local hod "Yes"	
		estadd local date "Yes"
		estadd local hours "All"
		estadd local sam "05/2015-08/2015"					
	est store C
	
*Column 4: Summer 2015 Diff + weather + date + HOD FE + pre-treatment use
	eststo: reghdfe euse ws elec_use_summer elec_use_annual elec_use_winter temp_60 temp_65 temp_70 temp_75 temp_80 temp_85 precip_1 if ym<=667, ///
			absorb(hour date) vce(cluster HHID)
	sum euse if e(sample) & ws==0
	scalar euse4 = r(mean)
		estadd scalar m_cont = scalar(euse4)
		estadd local weath "Yes"
		estadd local hod "Yes"	
		estadd local date "Yes"
		estadd local hours "All"
		estadd local sam "05/2015-08/2015"					
	est store D
		
*Column 5: Summer 2015 Diff weather + date + HOD FE - Peak Hours
	eststo: reghdfe euse ws temp_60 temp_65 temp_70 temp_75 temp_80 temp_85 precip_1 if ym<=667 & (hour>=15 & hour<=20), ///
			absorb(hour date) vce(cluster HHID)
	sum euse if e(sample) & ws==0
	scalar euse5 = r(mean)
		estadd scalar m_cont = scalar(euse5)
		estadd local weath "Yes"
		estadd local hod "Yes"	
		estadd local date "Yes"
		estadd local hours "Peak"
		estadd local sam "05/2015-08/2015"					
	est store E
	
*Column 6: Summer 2015 Diff + weather + date + HOD FE  + pre-treatment use- Peak Hours
	eststo: reghdfe euse ws elec_use_summer elec_use_annual elec_use_winter temp_60 temp_65 temp_70 temp_75 temp_80 temp_85 precip_1 if elec_use_summer!=. & ym<=667 & (hour>=15 & hour<=20), ///
			absorb(hour date) vce(cluster HHID)
	sum euse if e(sample) & ws==0
	scalar euse6 = r(mean)
		estadd scalar m_cont = scalar(euse6)
		estadd local weath "Yes"
		estadd local hod "Yes"	
		estadd local date "Yes"
		estadd local hours "Peak"
		estadd local sam "05/2015-08/2015"					
	est store F
	
esttab A B C D E F using $tables\eregs_$outputdate.tex, replace label ///
    booktabs b(a2) nonumber ///
	drop(temp_60 temp_65 temp_70 temp_75 temp_80 temp_85 precip_1) ///
    starlevels(* 0.10 ** 0.05 *** 0.01) ///
	title(Electricity Intent-to-Treat Effects ///
	     (Dependent Variable: Electricity Use (kWh/hr)) \label{eregs1})  ///
	cells(b(fmt(3) star) se(fmt(3) par)) ///
	note(Notes: The table reports intent-to-treat results from an OLS ///
		 regression of hourly electricity use on assignment to the treatment. ///
		 Columns 1 and 2 include all observations from May 15, 2015 to May ///
		 31, 2016. Columns 3 and 4 restrict the sample to the summer of 2015 ///
		 (May 15 to August 30). Columns 5 and 6 further limit the sample to ///
		 include only peak demand hours (3 PM to 8 PM). Standard errors are ///
		 clustered at the household. *, **, *** denote significance at the ///
		 10\%, 5\%, and 1\%level.) ///
	scalars("m_cont Mean Control Group Use" "weath Weather Controls" /// 
		    "hod Hour of Day FE" "date Calendar Date FE" ///
			 "hours Hours" "sam Sample") ///
	 mlabels((1) (2) (3) (4) (5) (6)) collabels(none) 

	 
********************************************************************************
*TABLE 2: PANEL A - WATER INTENT TO TREAT EFFECTS
********************************************************************************
eststo clear
*Column 1: Full Diff + weather + date + HOD FE
	eststo: reg wuse ws temp_60 temp_65 temp_70 temp_75 temp_80 temp_85 precip_1, vce(cluster HHID)
	sum wuse if e(sample) & ws==0
	scalar wuse1 = r(mean)
		estadd scalar m_cont = scalar(wuse1)
		estadd local weath "No"
		estadd local hod "No"	
		estadd local date "No"
		estadd local hours "All"
		estadd local sam "05/2015-05/2016"	
	est store A
	
*Column 2: Full Diff + weather + date + HOD FE + pre-treatment use
	eststo: reghdfe wuse ws w_use_summer w_use_annual w_use_winter temp_60 temp_65 temp_70 temp_75 temp_80 temp_85 precip_1, ///
			absorb(hour date) vce(cluster HHID)
	sum wuse if e(sample) & ws==0
	scalar wuse2 = r(mean)
		estadd scalar m_cont = scalar(wuse2)
		estadd local weath "Yes"
		estadd local hod "Yes"	
		estadd local date "Yes"
		estadd local hours "All"
		estadd local sam "05/2015-05/2016"	
	est store B
		
*Column 3: Summer 2015 Diff  + weather + date + HOD FE
	eststo: reg wuse ws temp_60 temp_65 temp_70 temp_75 temp_80 temp_85 precip_1 if  ym<=667, vce(cluster HHID)
	sum wuse if e(sample) & ws==0
	scalar wuse3 = r(mean)
		estadd scalar m_cont = scalar(wuse3)
		estadd local weath "No"
		estadd local hod "No"	
		estadd local date "No"
		estadd local hours "All"
		estadd local sam "05/2015-08/2015"					
	est store C
	
*Column 4: Summer 2015 Diff + weather + date + HOD FE + pre-treatment use
	eststo: reghdfe wuse ws w_use_summer w_use_annual w_use_winter temp_60 temp_65 temp_70 temp_75 temp_80 temp_85 precip_1 if ym<=667, ///
			absorb(hour date) vce(cluster HHID)
	sum wuse if e(sample) & ws==0
	scalar wuse4 = r(mean)
		estadd scalar m_cont = scalar(wuse4)
		estadd local weath "Yes"
		estadd local hod "Yes"	
		estadd local date "Yes"
		estadd local hours "All"
		estadd local sam "05/2015-08/2015"					
	est store D

*Column 5: Summer 2015 Diff Peak Hours + weather + date + HOD FE
	eststo: reg wuse ws temp_60 temp_65 temp_70 temp_75 temp_80 temp_85 precip_1 if ym<=667 & (hour>=15 & hour<=20), vce(cluster HHID)
	sum wuse if e(sample) & ws==0
	scalar wuse3 = r(mean)
		estadd scalar m_cont = scalar(wuse3)
		estadd local weath "No"
		estadd local hod "No"	
		estadd local date "No"
		estadd local hours "Peak"
		estadd local sam "05/2015-08/2015"					
	est store E
	
*Column 6: Summer 2015 Diff Peak Hours + weather + date + HOD FE + pre-treatment use
	eststo: reghdfe wuse ws w_use_summer w_use_annual w_use_winter temp_60 temp_65 temp_70 temp_75 temp_80 temp_85 precip_1 if ym<=667 & (hour>=15 & hour<=20), ///
			absorb(hour date) vce(cluster HHID)
	sum wuse if e(sample) & ws==0
	scalar wuse4 = r(mean)
		estadd scalar m_cont = scalar(wuse4)
		estadd local weath "Yes"
		estadd local hod "Yes"	
		estadd local date "Yes"
		estadd local hours "Peak"
		estadd local sam "05/2015-08/2015"					
	est store F
	
esttab A B C D E F using $tables\wregs_$outputdate.tex, replace label ///
    booktabs b(a2) nonumber ///
	drop(_cons temp_60 temp_65 temp_70 temp_75 temp_80 temp_85 precip_1) ///
    starlevels(* 0.10 ** 0.05 *** 0.01) ///
	title(Water Intent-to-Treat Effects (Dependent Variable: ///
		  Water Use (gals/hr)) \label{wregs1})  ///
	cells(b(fmt(3) star) se(fmt(3) par)) ///
	note(Notes: The table reports intent-to-treat results from an OLS ///
	     regression of hourly water use on assignment to the treatment. ///
		Columns 1 and 2 include all observations from May 15, 2015 to May ///
		31, 2016. Columns 3 and 4 restrict the sample to the summer of 2015 ///
		(May 15 to August 30). Columns 5 and 6 further limit the sample to ///
		include only peak demand hours (3 PM to 8 PM). Standard errors are ///
		clustered at the household. *, **, *** denote significance at the ///
		10\%, 5\%, and 1\%level.) ///
	scalars("m_cont Mean Control Group Use" "weath Weather Controls" /// 
		    "hod Hour of Day FE" "date Calendar Date FE" ///
			"hours Hours" "sam Sample") ///
	 mlabels((1) (2) (3) (4) (5) (6)) collabels(none) 		
	

********************************************************************************
*TABLE 3: TREATMENT EFFECT HETEROGENEITY AND OUTDOOR TEMPERATURE
********************************************************************************
*Generating watering days indicator
* Notes: - Before July 2014: 
*				No Watering Restrictions
*		 - Starting July 22, 2014: 
*				Tues/Thurs/Sat from 6/22/2014-10/31/2014 
*				Sat. from 11/1/2014-3/31/2015
*				Tues/Thurs/Sat from 4/1/2015-5/31/2015
*		 - Starting June 1, 2015
*				Tues/Sat from 6/1/2015-10/31/2015
*				Sat from 11/1/2015-3/31/2016
*				Tues/Sat from 4/1/2016-6/20/2016
*		 - Starting June 21, 2016
*				Tues/Thurs/Sat from 6/21/2016-10/31/2016
gen dow=dow(date) 
gen wat_day=0
	replace wat_day=1 if dow==2|dow==4|dow==6 & date>20169 & date<=20239
	replace wat_day=1 if dow==2|dow==6 & date>=20240 & date<=20392
	replace wat_day=1 if dow==6 & date>=20393 & date<=20544
	replace wat_day=1 if dow==2|dow==6 & date>=20545 & date <=20625
	replace wat_day=1 if dow==2|dow==4|dow==6 & date>=20626

*Merging historical water data	
merge m:1 HHID using $clean\BaselineWater_20171018.dta
	drop if _merge==2
	drop _merge
drop control WS HWS HWSPlus

*Ordering variables
order HHID date_h euse wuse ws temp_ws_* temp_6* temp_7* temp_8* 

eststo clear		
	
*Column 1: Electricity - Summer 2015 + Precip + Date FE + No Watering
eststo: reghdfe euse temp_65-temp_85 temp_ws_60-temp_ws_85 elec_use_summer elec_use_annual elec_use_winter precip_1 /// 
	if ym<=667 & elec_use_summer!=. & wat_day==0, a(date) vce(cluster HHID)
	test [temp_ws_65]=[temp_ws_85]
		estadd local dep "Electricity"	
		estadd local date "Yes"
		estadd local water "No Watering"		
		estadd local sam "5/15-8/15"					
	est store A
	
*Column 2: Electricity - Summer 2015 + Precip + Date FE + Watering
eststo: reghdfe euse temp_65-temp_85 temp_ws_60-temp_ws_85 elec_use_summer elec_use_annual elec_use_winter precip_1 /// 
	if ym<=667 & elec_use_summer!=. & wat_day==1, a(date) vce(cluster HHID)
	test [temp_ws_65]=[temp_ws_85]
		estadd local dep "Electricity"	
		estadd local date "Yes"
		estadd local water "Watering"		
		estadd local sam "5/15-8/31"					
	est store B
	
*Column 3: Water - Summer 2015 + Precip + Date FE + No Watering
eststo: reghdfe wuse temp_65-temp_85 temp_ws_60-temp_ws_85 w_use_summer w_use_annual w_use_winter precip_1 ///
		if ym<=667 & elec_use_summer!=. & wat_day==0, a(date) vce(cluster HHID)
	test [temp_ws_60]=[temp_ws_85]
		estadd local dep "Water"	
		estadd local date "Yes"
		estadd local water "No Watering"		
		estadd local sam "5/15-8/31"					
	est store C
	
*Column 4: Water - Summer 2015 + Precip + Date FE + Watering
eststo: reghdfe wuse temp_65-temp_85 temp_ws_60-temp_ws_85 w_use_summer w_use_annual w_use_winter precip_1 ///
		if ym<=667 & elec_use_summer!=.  & wat_day==1, a(date) vce(cluster HHID)
	test [temp_ws_60]=[temp_ws_85]
		estadd local dep "Water"	
		estadd local date "Yes"
		estadd local water "Watering"		
		estadd local sam "5/15-8/31"					
	est store D
		
esttab A B C D using $tables\regs_temp_$outputdate.tex, replace label ///
    booktabs b(a2) nonumber ///
    starlevels(* 0.10 ** 0.05 *** 0.01) ///
	title(Heterogeneity by Outdoor Temperature ///
	     (Dependent Variable: Electricity and Water Use) \label{regs_temp})  ///
	cells(b(fmt(3) star) se(fmt(3) par)) ///
	note(Notes: The table reports intent-to-treat effects disaggregated ///
		 by hours in which outdoor temperatures lie in each bin. Columns ///
		 1 to 3 report the results for electricity use and columns 4 to ///
		 6 report the results for water use. Day of Week indicates the ///
		 days included in our sample.  `All' includes all data for Summer ///
		 2015, `No Watering' restricts the sample to days of the week when ///
		 outdoor watering was banned, and `Watering' restricts the sample ///
		 to days of the week when outdoor watering was allowed.  All ///
		 regressions include controls for the temperature bins 65F-70F, ///
		 70F-75F, 75F-80F, 80F-85F and $\geq$ 85F, hourly outdoor ///
		 precipitation, and calendar date fixed effects. Standard errors ///
		 in parentheses are clustered at the household. *, **, *** denote ///
		 significance at the 10\%, 5\%, and 1\% level.) ///
	scalars("dep Dependent Variable" "date Calendar Date FE"  ///
			 "water Day of Week" "sam Sample") ///
	 mlabels((1) (2) (3) (4)) collabels(none) 	  
	  	
	
********************************************************************************
*FIGURE 3: ELECTRICITY INTENT TO TREAT EFFECTS OVER TIME
********************************************************************************
*Year-Month Indicators and Year-Month Interactions
tab ym, gen(ym_) 

 forval i = 1/13 {
	local j=663+`i'
	gen ym_ws_`i'=treat*ym_`i'
	label var ym_ws_`i' "WS, YM=`j'"
	disp `i'
 }
 *
compress

*Regression
reghdfe euse ym_ws_1-ym_ws_13 ///
	  elec_use_summer elec_use_annual elec_use_winter ///
	  temp_60 temp_65 temp_70 temp_75 temp_80 precip_1, ///
	  absorb(hour date)  cluster(HHID)
			
*Saving results
gen ws_ym=.
	label var ws_ym "Month"
gen euse_ym=.
	gen euse_ym_ub=.
	gen euse_ym_lb=.
forval i = 1/13 {
	local j=663+`i'
	replace ws_ym=`j' in `i'
	replace euse_ym=_b[ym_ws_`i'] in `i'
	replace euse_ym_ub=(_b[ym_ws_`i']+1.96*_se[ym_ws_`i']) in `i'
	replace euse_ym_lb=(_b[ym_ws_`i']-1.96*_se[ym_ws_`i']) in `i'	
}
*		
*Regression
reghdfe euse ym_ws_1-ym_ws_13 ///
	  elec_use_summer elec_use_annual elec_use_winter ///
	  temp_60 temp_65 temp_70 temp_75 temp_80 precip_1 if (hour>=15 & hour<=20), ///
	  absorb(hour date)  cluster(HHID)
gen euse_ym_peak=.
	gen euse_ym_peak_ub=.
	gen euse_ym_peak_lb=.
forval i = 1/13 {
	replace euse_ym_peak=_b[ym_ws_`i'] in `i'
	replace euse_ym_peak_ub=(_b[ym_ws_`i']+1.96*_se[ym_ws_`i']) in `i'
	replace euse_ym_peak_lb=(_b[ym_ws_`i']-1.96*_se[ym_ws_`i']) in `i'	
}
*	
preserve
keep ws_ym-euse_ym_peak_lb
drop if missing(ws_ym)

*Graph - TEs Over Time
format ws_ym %tm
twoway (line euse_ym_peak_ub ws_ym, lp(dash) lc(erose)) ///
		(line euse_ym_peak_lb ws_ym, lp(dash) lc(erose)) ///
       (connected euse_ym ws_ym, lc(edkblue) m(circle) mc(edkblue) msize(medium)  mfc(white) ) ///
	   (connected euse_ym_peak ws_ym, lc(cranberry) m(diamond) mc(cranberry) msize(medium)  mfc(white) ///
	   graphregion(color(white)) bgcolor(white)  ///
	   title("") xtitle("") ///  
	   ytit("Electricity ITT (kWh/hour)") ///
	   ylabel(-0.1(0.05)0.1, nogrid angle(hor)) ///
	   yline(0, lc(black) lp(solid)) ///
	   xlabel(664(4)676, angle(45)) ///
	   legend(order(3 "All Hours" 4 "Peak Hours") size(small) ///
	   region(lcolor(white)) cols(1) ring(0) position(5))  )
graph export $figs\euse_ym_controls_$outputdate.png, width(4000) replace	
restore
drop ws_ym-euse_ym_peak_lb
drop ym_* _est*

********************************************************************************
*FIGURE 4A: HETEROGENEITY BY HOUR OF DAY - ELECTRICITY
********************************************************************************

*Creating Hourly Indcators and Interactions 
forvalues i=0(1)23{
gen hod_`i'=0
	replace hod_`i'=1 if hour==`i'
	
gen hodxws_`i'=ws*hod_`i'
}
*
order hod_* hodxws_*, last

gen ws_hour=.
	label var ws_hour "Hour"
gen elec_h=.
	gen elec_h_ub=.
	gen elec_h_lb=.
gen elec_base_h=.
compress

*Regression - Summer 2015 with Weather Controls and Date FEs
reghdfe euse hod_1-hod_23 hodxws_* elec_use_summer elec_use_annual elec_use_winter temp_60 temp_65 temp_70 temp_75 temp_80 temp_85 precip_1 if ym<=667, ///
	absorb(date) vce(cluster HHID)

forvalues i=0(1)23{
	local j=`i'+1
	sum euse if e(sample) & ws==0 & hour==`i'	
		replace ws_hour=`i' in `j'
		replace elec_h=_b[hodxws_`i'] in `j'
		replace elec_h_ub=(_b[hodxws_`i']+1.96*_se[hodxws_`i']) in `j'
		replace elec_h_lb=(_b[hodxws_`i']-1.96*_se[hodxws_`i']) in `j'
		replace elec_base_h=r(mean) in `j'
	display `i'
}
*
preserve
keep ws_hour-elec_base_h
drop if missing(ws_hour)

*Graph - Hourly TEs & Avg. Use
twoway rcap elec_h_ub elec_h_lb ws_hour, lstyle(ci) lc(erose) || ///
       scatter elec_h ws_hour, m(diamond) mc(cranberry) msize(medium) /// 
	   graphregion(color(white)) bgcolor(white) ///
	   title("") xtitle("Hour of Day") ///  
	   ytit("Electricity ITT (kWh/hour)") ///
	   legend(off) ///
	   xlabel(0(2)23) ///	   
	   yline(0, lcolor(black)) ///
	   ylabel(-0.1(0.05)0.01, nogrid)  ///
	   name(elec_1a, replace) fysize(70) 
twoway bar elec_base_h ws_hour, xlabel(0(2)23) ///
	   color(cranberry*0.6) graphregion(color(white)) ///
	   xtitle("Hour of Day") ytitle("Avg. Use") ///
	   bgcolor(white) ylabel(0(1)2, nogrid) ///
	   name(elec_1b, replace) fysize(30)
graph combine elec_1a elec_1b, cols(1)  iscale(1) /// 
              graphregion(color(white))  imargin(0 0 0 0) 
graph export $figs\euse_sum_hour_combine_$outputdate.png, width(4000) replace	
restore
drop ws_hour-elec_base_h


********************************************************************************
*FIGURE 4B: HETEROGENEITY BY HOUR OF DAY - WATER
********************************************************************************
gen ws_hour=.
	label var ws_hour "Hour"
gen wtr_h=.
	gen wtr_h_ub=.
	gen wtr_h_lb=.
gen wtr_base_h=.

*Regression - Summer 2015 with Weather Controls and Date FEs
reghdfe wuse hod_1-hod_23 hodxws_* w_use_summer w_use_annual w_use_winter temp_60 temp_65 temp_70 temp_75 temp_80 temp_85 precip_1 if ym<=667, ///
	absorb(date) vce(cluster HHID)

forvalues i=0(1)23{
	local j=`i'+1
	sum wuse if e(sample) & ws==0 & hour==`i'	
		replace ws_hour=`i' in `j'
		replace wtr_h=_b[hodxws_`i'] in `j'
		replace wtr_h_ub=(_b[hodxws_`i']+1.96*_se[hodxws_`i']) in `j'
		replace wtr_h_lb=(_b[hodxws_`i']-1.96*_se[hodxws_`i']) in `j'
		replace wtr_base_h=r(mean) in `j'
	display `i'
}
*
preserve
keep ws_hour wtr_h wtr_h_ub wtr_h_lb wtr_base_h
drop if missing(ws_hour)

*Graph - Hourly TEs & Avg. Use
twoway rcap wtr_h_ub wtr_h_lb ws_hour, lstyle(ci) lc(erose) || ///
       scatter wtr_h ws_hour, m(circle) mc(edkblue) msize(medium) /// 
	   graphregion(color(white)) bgcolor(white) ///
	   title("") xtitle("") ///  
	   ytit("Water ITT (gals/hr)") ///
	   legend(off) ///
	xlabel(none) xtick(0(2)23) ///	   
	yline(0, lcolor(black)) ///
	ylabel(-6(2)2, nogrid) ///
	name(wtr_1a, replace) fysize(70) 
twoway bar wtr_base_h ws_hour, xlabel(0(2)23) ///
	   color(edkblue*0.6) graphregion(color(white)) ///
	   xtitle("Hour of Day") ytitle("Avg. Use") ///
	   bgcolor(white) ylabel(0(15)30, nogrid) ///
	   name(wtr_1b, replace) fysize(30)
graph combine wtr_1a wtr_1b, cols(1)  iscale(1) /// 
              graphregion(color(white))  imargin(0 0 0 0) 
graph export $figs\wuse_sum_hour_combine_$outputdate.png, width(4000) replace
restore
frame change data_main
frame drop main_regs 